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Abstract. A knowledge of the recombination time on the grain surfaces 
has been a major obstacle in deciding the production rate of molecular hydro- 
gen and other molecules in the interstellar medium. We present a numerical 
study to compute this time for molecular hydrogen for various cloud and grain 
parameters. We also find the time dependence, particularly when a grain is 
freshly injected into the system. Apart from the fact that the recombination 
times seem to be functions of the grain parameters such as the activation barrier 
energy, temperature etc, our result also shows the dependence on the number 
of sites in the grain S and the effective accretion rate per site a s of atomic 
hydrogen. Simply put, the average time that a pair of atomic hydrogens will 
take to produce one molecular hydrogen depends on how heavily the grain is al- 
ready populated by atomic and molecular hydrogens and how fast the hopping 
and desorption times are. We show that if we write the average recombination 
time as T r ~ S a /An, where, Ah is the hopping rate, then a could be much 
greater than 1 for all astrophysically relevant accretion rates. Thus the average 
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formation rate of Hi is also dependent on the grain parameters, temperature 
and the accretion rate. We believe that our result will affect the overall rate of 
the formation of complex molecules such as methanol which require successive 
hydrogenation on the grain surfaces in the interstellar medium. 

Keywords : Molecular cloud - star formation - grain chemistry - numerical 
simulations 



1. Introduction 



It has long been suggested that the dust grains play a major role in the formation of 
molecular hydrogen in the interstellar medium (ISM) (Gould & Salpeter, 1963). Consid- 
erable studies were made since then to understand the real physical processes which are 
taking place both theoretically (e.g., Hollcnbach, Werner & Salpeter, 1971; Takahashi, 
Matsuda & Nagaoka, 1999; Biham et al. 2001) as well as experimentally (e.g., Pirroncllo 
ct al. 1997a,b, 1999). More recently, Biham et al. (2001), and Green et al. (2001) have 
computed H2 production rate by physisorption. It was found that a significant produc- 
tion is possible in cooler (<~ 10 — 25K) clouds. Cazaux & Tielens (2002, 2004) use both 
physisorption and chemisorption, to demonstrate that H2 production is possible at high 
temperatures (~ 200 — 400K) also. The goal is to study the rate at which the H atoms 
combine together on the surface of the grains to form H2 and then they are desorbed into 
the gas phase to react with other atoms. When compared with the average mass fractions 
of various molecular species obtained through gas phase reactions (see, Chakrabarti & 
Chakrabarti, 2000ab; Das et al. 2006), it was found that the observed abundances of 
more complex species, such as methanol, are much higher. It is possible that methanol 
as well as its precursors also have to be formed on grain surfaces through successive 
hydrogenation. Our finding for molecular hydrogen has thus important bearings on the 
formation of more complex molecules on grains. These molecules would, in turn, desorb 
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into the gas phase and would be expected to produce more complex species such as amino 
acids in due course. 

One of the most challenging problems is to determine the average rate at which 
the recombination of atomic hydrogen takes place on a grain surface. In theoretical 
investigations which are prevalent in the subject (See, Acharyya and Chakrabarti, 2005; 
hereafter Paper I; Acharyya, Chakrabarti and Chakrabarti, 2005; hereafter ACC05), the 
diffusion rate Ah (inverse of the diffusion time T<j = 1/A h) is divided by S, the number 
of sites on the grain surface (e.g., Biham et al. 2001) to get the recombination rate. The 
argument for reducing the rate by a factor of S is this: on an average, there are S 1 ^ 2 
number of sites in each direction of the grain. Since the hopping is random, it would take 
square of this, i.e., S number of hopping to reach a distance located at S 1 ^ 2 sites away, 
where, on an average, another H is available. Thus, the effective recombination rate was 
chosen to be Ah/S. It is an empirical factor and needs more careful treatment. In our 
present paper, we replace S by an 'unknown' quantity S' = S a (t), where a(t) may be 
time dependent (if the grain and cloud parameters change) and it could also deviate from 
unity. Let a — ao when a steady state is reached. Higher the accretion rate, lesser should 
be the value of ao as the effective surface area S' (t — > oo) gets smaller and smaller. The 
opposite is true for smaller accretion rates. In fact, in the limit, if the accretion rate 
is so low, that a lone H sweeps around the grain several times to find another H, one 
would get ct{) > 1 since the effective site number is higher than S. Using our simulation, 
we determine how the effective site number deviates from S, one way or the other when 
the accretion rate is varied. Our result is likely to have important consequences for the 
formation of other hydrogenated species, such as water, methanol on grain surfaces. This 
will be discussed elsewhere. Some preliminary results with steady state ao have been 
presented in Das et al. (2005) and Chakrabarti et al. (2006). In the current paper, 
we discuss the time and temperature dependence of a(t) and studied the cases for more 
varied astrophysically important accretion rates. 

In the next Section, we present the modified typical equations which govern the 
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molecular hydrogen production rate on a grain surface. Incorporating the physical aspects 
of these equations, we perform a numerical simulation to determine the numbers on the 
grains. In §2, we present the procedure for the simulation and in §3 we present the results 
for two types of commonly used grains, namely olivine and amorphous carbon. We show 
how a(t) depends on time and how it settles into a number (generally > 1) when steady 
state is reached. Finally in §4, we present our concluding remarks. 



2. Procedure of numerical simulation 



The relevant equations which are generally solved on grain surfaces have been presented 
in Paper I (Eqs. 2a-d) and we do not repeat them here. However, for Monte-Carlo 
simulation we need to modify these equations. Since in our simulation we expect to get 
the rate of diffusion (T r (t) = S/Ah) of H exactly, we can assume that the recombination 
time at any instant t could be written in the form: 

T r (t) = S'/A H , (1) 

where, S' is an 'effective surface area' which we may be written as S' = S a (t). Here, 
(f) H = F H (1 — f grh — fg r h2)- F H is the accretion rate of H and F H = A < v > N H , A 
being the area of a grain, < v > is the average velocity and N H is the number density of 
H in the gas phase. Our goal would be to compute the steady state value ceo = a (t — * oo) 
for various grains (type and size) at various accretion rates and grain temperatures and 
to check if a ~ 1. Thus, instead of Eqs. 2a of Paper I, the actual governing equation for 
atomic hydrogen should be, 

^f- = ^ H - W H n H - 2(A H /S a )n 2 H , (2a) 
at 

along with equation governing the rate of production of H 2 given by, 

^ = »(A H /S a )n% - W H2 n H2 , (2b) 

where, Ah is the hopping rate given by ;/exp(— Eo/kbT), Wr is the desorption co-efficient 
of hydrogen atom given by, j/exp(— Ei/kbT), Wh2 is the desorption co-efficient of hydro- 
gen molecule given by, v exp(— E 2 / kbT) . Here, Eq is the activation barrier energy for 
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diffusion of H atom, E\ is the activation barrier energy for desorption of H atom and 
E2 is the activation barrier energy for desorption of H2 molecule, v is the vibrational 
frequency given by, be the number of H atoms, and Uh2 be the number 

of H2 on a grain at time t, s ~ 10 14 ~ 16 is the surface density of sites on a grain, ruin is 
the mass of the H atom and Ed is the binding energy The parameter /i represents the 
fraction of H2 molecules that remains on the surface upon formation while (1-/U.) fraction 
is desorbed due to the energy released in the recombination process. The H2 production 
rate Rh2 in the gas due to grain is then given by, 

R H 2 = (1- v)(A H /S a )n 2 H + W H 2n H 2. (3) 

The values of the activation barrier energies E , E\ and E 2 are taken from Katz et al. 
(1999) and are given by, E = 24.7 mcV, E x = 32.1 mcV and E 2 = 27.1 mcV for olivine 
and E = 44 mcV, E\ = 56.7 meV and E 2 = 46.7 mcV for amorphous carbon grains. 
For olivine, [i = 0.33 and for amorphous carbon /i = 0.413 was used. 

For the sake of simplicity, we assume each grain to be square in shape having S = n 2 
number of sites (square lattice). We use periodic boundary condition so that the lattice 
behaves like a closed surface. We anticipate that the result would be insensitive to the 
actual nature of the lattice when simulation is done for a long time. We use monolayer on 
the grain surface. We choose the minimum time step to be th — 1/Au, the hopping time 
and advance the global time by this step. If 4>h is the effective accretion rate (as defined 
above), after every time step, i.e., after every 1/Ah seconds, 4>h/Ah number of hydrogen 
atoms are dropped on the grain. If <f>u is too low so that 4>h/Ah < 1, then, clearly, one 
H is dropped after every Ah/4>h steps. The exact site at which one atom is dropped 
is obtained by a pair of random numbers (R Xl R yi R Xl R y < 1) obtained by a random 
number generator. This pair would place the incoming hydrogen at (z,j)th grid, where, 
i and j are the nearest integers obtained using the Int function: i = int(R x * n + 0.5) 
and j = int(R y * n + 0.5). Each atom starts hopping with equal probability in all four 
directions. This was decided by another random number. When, during the hopping 
process, one atom enters a site which is already occupied by another atom, we assume 
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that a molecule has been formed and increase the number of H 2 by unity. However, when 
the atom enters a site occupied by an H 2l another random number is generated to decide 
which one of the other neatest sites it is going to occupy. 

Thermal evaporation of H and H 2 from a grain surface are handled in the following 
way. Since Wr is the desorption rate for H, one atom is supposed to be released to the 
gas phase after every 1/Wh seconds. We generate a random number R t for every H 
present on the grain and release (at each time step, i.e., 1/A h seconds) only those for 
which R t < Wh/Ah- Similar procedure is followed for the evaporation of H 2 for which 
the criterion for evaporation was R t < Wh2/Ah- 

Due to spontaneous desorption, a factor of of riH2 is lost to the gas phase. Here 
too, a random number R s is generated for each newly formed (within that time step) H 2 
present on the grain. Those which satisfy R s < (1 — 11) are removed to the gas phase. 

We continue our simulation for more than 10 9 s. After an initial transient period, n# 
and UH2 reaches a quasi-steady state with some fluctuations. In a steady state, ao can 
be calculated using Eq. 2a, i.e., 

9 Ann 2 

ao = log(- y^)/log(S). (4) 

<pH ~ W H n H 

In this context, one important parameter [3 may be defined as the 'catalytic capacity' of 
a grain (Chakrabarti et al. 2006) which measures the efficiency of the formation of H 2 
on that grain surface for a given pair of H residing on it. Let 8Nh 2 be the number of H 2 
formed in St time. Since two hydrogen atoms are required to create one H 2 , the average 
rate of creation of one H 2 per pair of H atom would be given by, 

We identify the inverse of this rate with the average formation rate given by, 

7>(t) = S^/A H . (56) 

Thus, 

S f3 (t) = A H / <A H1 >. (5c) 
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This yields /3(f) as a function of time as, 



(3 = log(A H / < A H1 >)/log(S). 



(6) 



Here too we can define /3 = /3(f — > oo). Our goal would be to see if flo thus obtained 
actually varies with grain parameters, accretion rate and temperatures. Intuitively, it 
should. This is because, when the accretion rate is very low, finding a second H on 
the grain would be difficult and it can take several sweeps of the grain surface. Thus 
/?o > 1 would be a possibility. On the other hand, when the rate is very high, many sites 
would be occupied and another H would be met along any direction and /3o could come 
down to 1/2 or even less. Indeed, this is what we see as well. We have carried out the 
simulation with and without spontaneous desorption taken into account. In presence of 
the spontaneous desorption, we remove Hi by generating a random number as soon as 
one Hi is formed and checking if it is less or more compared to /i. If less, Hi remains on 
the grain, else it is taken out to the gas phase. The detailed results are discussed in the 
next Section. 



Fig. 1 gives the variation of nn and uh2 residing on the grain surface since the beginning 
of the simulation. We chose a grain with 10 4 sites at 8K. For clarity, we plot average 
numbers in every ~ 5.5 x 10 6 s bin after the initial transient period of ~ 2 x 10 6 s is over. 
The simulation was carried out near about 10 10 s. The effective accretion rate is assumed 
to be <j)R = 7.98 x 10 _4 s _1 . In this case, the steady state has clearly been reached by 
this time and we can compute a by using Eq. 4. This is done below. 

3.1 Olivine grains 

We start with the olivine grains which are kept at 8 K. The evaporation rate of Hi is 
Wh 2 = 0.834844 x 10 -5 s^ 1 which corresponds to a time scale of <~ 119782.8s and similarly 
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Figure 1. : Variation of the number of H and Hi on an olivine grain of 10 4 sites kept at 8K 
and exposed to an effective accretion rate per site of H 7.98 x 1CF 8 s _1 . We carried out our 
simulation near about 10 10 s. For cleanliness, after the initial transient period of about 2 x 10 6 s 
is over, we took time average at every 5.5 x 10 6 s before plotting the numbers. 

the evaporation rate of H is Wh = 0.59 x 10 _8 s _1 which corresponds to ~ 169585908s. 
The hopping time for hydrogen on this grain is 3680.58s (data taken from Katz et al. 
1999). 

In Fig. 2a, we present the computed a® (Eq. 4) as a function of a s - the effective 
accretion rate per site. The solid, dot-dashed and the dashed curves are for 5* = 10 4 , 
9 x 10 4 and for 10 6 sites respectively. No spontaneous desorption has been included (i.e., 
/i = 1). We extrapolated our curve to very low accretion rates which would have taken a 
very long computation time, just to show the trend of the result and extreme conditions. 
We note that «o is generally higher than unity in the region of our interest. This is 
because Hs are scattered few and far between and it takes a longer time (generally more 
than one sweeping) for one H to locate another. a$ monotonically drops as the accretion 
rate goes up. In Fig. 2b, we show the variation of (3q (Eq. 6) for the same case. Here too, 
we see that [3o is very high compared to unity for low rates, but becomes ~ 0.5 or lower 
for higher rates as expected. Note that «o and (3q go down with increasing site number S 
also. Since for a smaller grain, the possibility of getting it filled at a high rate is higher, 
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Figure 2. (a-b): Variation of (a) Qo and (b) /3o as a function of a s , the effective accretion rate 
per site for various olivine grains kept at 8K. The dashed, dot-dashed and solid curves are for 
S — 10 6 , 9 x 10 4 , and 10 4 respectively. «o is clearly a function of the accretion rate. For rates 
relevant in molecular clouds ao and f3o are much larger than unity. For very high rates (3q comes 
down to 0.5 or lower. The deviation is highlighted using dotted curves by extrapolating at very 
low accretion rates. 

one would have expected an opposite result. However, it is to be remembered that for a 
larger grain, the accretion rate itself (4>h — Sa s ) is also large. Hence the plots are to be 
compared carefully. For instance, the result of a s = 0.0005 x 10 _8 s _1 for 9 x 10 4 sites 
is to be compared with that of a s ~ 0.0045 x 10 _8 s _1 for 10 4 sites in order to make a 
meaningful comparison. In any case, for reasonable 4>h values with number densities up 
to 10 6 cm~ 3 , the relevant a s would be below 10" 6 where ao > 1 in general. 
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Figure 3. (a-d): In (a-b) snapshots of the grain surface with H (hollow squares) and H2 (filled 
squares) at two arbitrarily chosen times (a) 8 x 10 8 s and at (b)10 9 s. Here an olivine grain (at 
8K) with 900 sites has been chosen. This is bombarded with an accretion rate per site of H 
3.02 x 10~ 7 per sec. No spontaneous desorption has been assumed here. In (c-d) spontaneous 
desorption has been included and plotted for the same time as before. Thus, numbers of H2 
residing on the grain at any instant are lesser. 

In Fig. 3(a-b), we present snapshots of the occupancy of H (hollow squares) and 
Hi (filled squares) at two instants of time on the grain containing only 900 sites at two 
arbitrarily chosen times (a) t=8 x 10 8 s (b) t=10 9 s respectively. Here an olivine grain at 
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Figure 4. : Temperature dependence of cto for the olivine grains at 10K (solid), 9K (dot-dashed) 
and 8K (dashed). The deviation is highlighted using dotted curves by extrapolating at very low 
accretion rates. 

(8K) with 900 sites has been chosen. This grain is facing an accretion rate of Fh / S of 
3.02 x 10~ 7 s _1 per site. No spontaneous desorption has been assumed here. This is to 
be compared with Fig. 3(c-d) for the same simulation and for the same time as before 
when the spontaneous desorption has been included. The numbers of Hi are fewer since 
after formation of H2 molecules some part of the H2 are spontaneously desorbed in to 
the gas phase. 

When temperatures of the grain is increased, all the rates go down exponentially. As 
a result, we expect ao to rise with temperature for a given accretion rate. In Fig. 4 we 
show this behaviour for T = 8, 9, 10K respectively for olivine grains. We choose S = 10 4 
in this case. This behaviour affects the recombination efficiency r\ as defined by, 

n A (7) 

itself. In Fig. 5, we compare the temperature dependence of r\ as obtained from the rate 
equations (solid) curve with that obtained from our simulations. The accretion rate per 
site of 1.8 x 10~ 9 scc^ 1 was used (same in both the cases). We note that for T < 7.5K 
the simulation results are higher and for T ^ 7.5K the simulation results are lower. This 
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Figure 5. : A comparison between the recombination efficiency obtained from the rate equation 
(solid) and that obtained from our simulation (dashed). We use the accretion rate per site 
1.8 x 1CP 9 per second for a grain of diameter 0.1/xm. The difference can be attributed to the 
temperature dependence of the ao as shown in Fig. 4 above. 
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Figure 6. : A comparison of the simulation results (dark circles) with those obtained from 
analytical considerations (dashed curves) when suitable modification of the average recombina- 
tion rate is made. An Olivine grain of 10 4 sites at a temperature of 8K has been chosen in this 
comparison. Dotted curves are drawn using analytical results for o?o extrapolated to very low 
accretion rates. 
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is because ao itself is strongly temperature dependent as shown in Fig. 4 while the rate 
equation uses S = S' (i.e., ao = 1) for all temperatures. 

In Fig. 6 we compare our results with the analytically obtained results from the 
rate equation method provided correct S' was chosen (Eqs. 2ab). The simulation results 
are shown by the dark circles and those obtained from the analytical considerations are 
shown by the dashed curves. An Olivine grain of 10 4 sites at a temperature of 8K has 
been chosen in this comparison. Dotted curves are drawn using analytical results for a 
extrapolated to very low accretion rates. 

It is interesting to compare the results of our simulation with those obtained from 
the analytical considerations with and without our ao factor. In Table 1, we present this 
comparison. We take an Olivine grain of 10 4 sites at 8K and vary the accretion rates. 
In Column 1, we give the accretion rate per site of the grain. In Column 2 we present 
the coefficient ao which we derive from our simulation. In Columns 3-5, we present the 
number of H as obtained by our simulation and the modified equation (Eq. 2a) and 
the standard equations (Paper I) respectively. Columns 6-8, we present similar results 
for H 2 - We find that our simulation matches more accurately with the analytical results 
provided 5" is chosen as the surface area. If the standard equation is used, the deviation 
is very significant. Indeed, the number of H on the grain could be roughly half as much 
when simplistic analytical model is used. What observe is that on the grains we tend to 
have more H and less H2 than what analytical work suggests. 

3.2 Amorphous carbon grains 

Carbon grains produce significant H2 at a higher temperature than olivine grains because 
all the barrier energies are higher. We plot in Fig. 7a the variation of ao with a s at 
temperature 14K. The nature of variation of ao remains the same, namely, ao goes down 
with a s . The solid, dot-dashed and the dashed curves are for £ = 10 6 , 9 x 10 4 and for 
10 4 sites respectively. In Fig. 7b, we show the variation of (3q and as expected its value 
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Table 1. Comparison of H and H2 abundances in various methods 
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simulation 
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ceo = 1 


simulation 


CtQ 7^ 1 




ceo = 1 




6.79 x 1CT 7 


1 


04 


403.11 


407.80 


340.57 


361.71 


383.82 




377.48 




2.72 x 1CT 7 


1 


05 


275.78 


278.02 


219.31 


149.96 


158.16 




156.53 




2.72 x 10~ 8 


1 


07 


99.03 


99.31 


70.35 


15.40 


16.14 




16.11 




1.36 x 1CT 8 


1 
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72.27 


72.42 


49.80 


7.71 


8.08 




8.07 




5.43 x 1CT 9 


1 
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47.37 


31.52 


3.11 


3.23 




3.23 




2.72 x 10~ 10 


1 


11 


11.62 


11.62 


7.01 


1.49 x 10- 1 


1.59 x 10 


-1 


1.60 x 10 


-1 
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7.46 
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5.81 x 10~ 2 


6.24 x 10 


-2 


6.35 x 10 


-2 


2.72 x 10~ n 


1 


13 


3.86 
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2.18 


1.31 x 10~ 2 


1.49 x 10 


-2 


1.55 x 10 


-2 


*2.04 x 10~ n 
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1.11 x 10 


-2 


1.17 x 10 


-2 


*2.92 x 10~ 12 


1 
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1.32 


1.30 


0.70 


1.3 x 10~ 3 


1.5 x 10- 


3 


1.6 x 10" 


3 


*4.16 x 10~ 13 


1 


16 


0.52 


0.51 


0.26 


2.00 x 10~ 4 


2.00 x 10 


-4 


2.00 x 10 


-4 


*5.95 x 10~ 14 


1 


17 


0.20 


0.20 


0.10 


2.11 x 10- 5 


2.81 x 10 


-5 


3.13 x 10 


-5 



* represents the extrapolated value. 

can become as low as 0.5 for very large accretion rate. These results are representative 
as they are strongly temperature dependent as in the case of olivine (see, Fig. 4). 



4. Concluding remarks 

In this paper, we studied the behaviour of the average recombination time of atomic 
hydrogens on grain surfaces as a function of time, temperature, accretion rate and grain 
parameters. In the literature, it is simply assumed that the recombination time is grain 
site number divided by the diffusion rate. However, we find that this simplistic assumption 
is not valid especially when the accretion rate of H on the grain is very low or very high. 
For very low rates, an atomic hydrogen may have to sweep a grain several times before 
meeting another atom to form a molecule. For very high rates, the grain surfaces could be 
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Figure 7. (a-b): Variation of (a) ao and (b) /3o as a function of a 3 , the effective accretion rate 
per site, for various amorphous carbon grains kept at 14K. The solid, dot-dashed and dashed 
curves are for S = 10 4 , 9 x 10 4 , 10 6 respectively The exponents ao and flo are strong functions 
of the accretion rate. For low rates f3o is higher than unity, while for higher rates it is close to 
0.5 or even lower. The deviation is highlighted using dotted curves by extrapolating at very low 
accretion rates. 

partially filled with H and H 2 and search for another H need not take S number of hops. 
We found that on an average, ao, the index which determines the average recombination 
time, could be grater than 1 for very low rates. We also find that ao depends very 
strongly on the nature of the grains as well as the temperature and grain site numbers. 
We show how the number of H and H2 change with time and determine the time scales in 
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which quasi-cquilibrium is reached. Thus the recombination time is a complex function 
of these parameters. We defined another index (3 (Chakrabarti et al. 2006) which is a 
measure of the average formation rate of H^. This index also is also a strong function of 
accretion rate per site. In fact, for very high rates the index may go down to 0.5 or even 
lower. A comparison of the recombination efficiency as obtained from our procedure with 
that obtained from the rate equation clearly shows a deviation easily attributable to the 
temperature dependence of a in a realistic case. 

The results we obtained are likely to be astrophysically important because of several 
factors. The physical reason of the dependence of the average recombination time is so 
generic that one needs to revise all the abundances made on the grain surfaces. We find 
that even for H, we have more H left over on grains, while less H 2 is produced and 
released. Thus, analytical considerations will always over produce molecules obtained 
on grain surfaces. Similarly, if, for example, methanols are produced through successive 
hydrogenation process on a grain surface, at each step such corrections would be required 
and the final abundance would be greatly affected. This will be done in near future. 
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